home *** CD-ROM | disk | FTP | other *** search
- //-------------------------------------------------------------------//
-
- // Synopsis: Polar decomposition.
-
- // Syntax: P = poldec ( A )
-
- // Description:
-
- // Poldec computes a matrix U of the same dimension as A, and a
- // Hermitian positive semi-definite matrix H, such that A =
- // U*H. U has orthonormal columns if m>=n, and orthonormal rows
- // if m<=n. U and H are computed via an SVD of A. U is a nearest
- // unitary matrix to A in both the 2-norm and the Frobenius norm.
-
- // Poldec returns a list with elements `u' and `h'.
-
- // Reference:
- // N.J. Higham, Computing the polar decomposition---with applications,
- // SIAM J. Sci. Stat. Comput., 7(4):1160--1174, 1986.
-
- // (The name `polar' is reserved for a graphics routine.)
-
- // This file is a translation of poldec.m from version 2.0 of
- // "The Test Matrix Toolbox for Matlab", described in Numerical
- // Analysis Report No. 237, December 1993, by N. J. Higham.
-
- //-------------------------------------------------------------------//
-
- poldec = function ( A )
- {
- local (A)
-
- m = A.nr; n = A.nc;
-
- svda = svd(A); // Economy size.
- svda.sigma = diag(svda.sigma);
-
- U = svda.u*svda.vt;
- H = svda.vt'*svda.sigma*svda.vt;
- H = (H + H')/2; // Force Hermitian by taking nearest Hermitian matrix.
-
- return << u = U ; h = H >>;
- };
-